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ABSTRACT 

CNl ■ We use an 850-/im SCUBA map of the HDF to study the dust properties of optically- 

selected starburst galaxies at high redshift. The optical/IR data in the HDF allow 
a photometric redshift to be estimated for each galaxy, together with an estimate of 
the visible star-formation rate. The 850-/im flux density of each source provides the 
complementary information: the amount of hidden, dust-enshrouded star formation 
activity. Although the 850-/im map does not allow detection of the majority of indi- 
| vidual sources, we show that the galaxies with the highest UV star-formation rates 

, are detected statistically, with a flux density of about Ssso = 0.2 mJy for an apparent 

UV star- formation rate of 1 ft. _2 M©yr _1 . This level of submillimetre output indicates 
that the total star-forming activity is on average a factor of approximately 6 times 
larger than the rate inferred from the UV output of these galaxies. The general popu- 
lation of optical starbursts is then predicted to contribute at least 25% of the 850- um 
^—j | background. We carry out a power-spectrum analysis of the map, which yields some 

evidence for angular clustering of the background source population, but at a level 
lower than that seen in Lyman-break galaxies. Together with other lines of argument, 
particularly from the NICMOS HDF data, this suggests that the 850-//m background 
originates over an extremely wide range of redshifts - perhaps 1 < z < 6. 
K*N ' 
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1 INTRODUCTION 

The observational ingredients for an understanding of galaxy 
formation are arguably now largely in place. In particu- 
lar, the Lyman-break technique has allowed the selection 
of galaxies over the redshift range 2 < z < 5 (Steidel et 
al. 1996). The combination of ground-based and HST data 
has allowed such samples to probe a wide range of opti- 
cal/UV luminosities, so that the statistical properties of the 
high-redshift galaxy population are increasingly well-known. 
This knowledge is conveniently summarized in the form of 
the global star- formation rate as a function of redshift. This 
function rises steeply between z = and 2 = 1 (Lilly et al. 
1996); initial indications from the HDF were that there was 
then a maximum, followed by a decline at z > 2 (Madau 
et al. 1996). More recent Lyman-break data, however, have 
cast doubt on this claim: it now seems more likely that the 
UV luminosity density remains approximately constant for 

2 < z < 5 (Steidel et al. 1999). Nevertheless, since there 



is little cosmological time at z > 5, these results are of- 
ten taken to suggest that we have seen almost all the star- 
formation that ever occurred. Provided we can understand 
the masses of the galaxies that are involved (as is starting 
to be possible though Balmer-line spectroscopy; Pettini et 
al. 1998), the entire history of the galaxy- formation process 
is potentially open to view. 

The main objection to this optimistic picture is well 
known: star formation in the local universe is often accom- 
panied by dust. There is thus the possibility that a large, 
or even dominant, component of the star-forming history 
is not revealed by optical/UV data. It has long been clear 
that this question would only be settled by studies in the 
far-IR or sub-mm bands, and great progress is now being 
made in this area. The total amount of energy released from 
young stars and reprocessed into emission from cool dust is 
constrained by the detection of the background radiation at 
A > 150 /im (Puget et al. 1996; Schlegel, Finkbeiner & Davis 
1998; Fixsen et al. 1998; Hauser et al. 1998). The resolution 
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of this background into point sources has been made pos- 
sible with the SCUBA imager on the James Clerk Maxwell 
Telescope (Holland et al. 1999). A number of groups have 
used this facility in order to map the sub-mm sky down 
to an rms sensitivity below 1 mjy at 850 /jm (Smail et al. 
1997; Hughes et al. 1998; Blain et al. 1999; Barger, Cowie 
& Sanders 1999; Eales et al. 1999). 

The observations of Hughes et al. (hereafter H98) are 
the deepest of existing datasets, and also have the advan- 
tage that they cover the area of the northern Hubble Deep 
Field (Williams et al. 1996). The galaxy population in this 
region of the sky has been studied in great detail, and so 
there is the potential to extract considerable information 
from the sub-mm data. The initial study by H98 considered 
only possible optical counterparts for the brightest 5 850- /im 
sources. Here, the intention is to dig much more deeply into 
the background population. The optical HDF data tell us 
the location of a large number of starburst galaxies, and it 
should be possible to measure in a statistical way their sub- 
mm emission. In this way, we will learn about the fraction 
of the star-forming activity in these galaxies that is hidden 
from view. The existence of photometric redshifts means 
that we can do this as a function of redshift; this measure- 
ment is critical, because it will tell us whether or not to 
believe the indications from the optical/UV data that most 
of the star-formation in the universe has been seen. 

In fact, there are already indications that a large frac- 
tion of the star-forming activity in the universe is not visible 
in the UV. From the first general detection of luminous in- 
frared galaxies (e.g. Soifer et al. 1987), it has increasingly 
become clear that the majority of the star-forming activity 
at low redshift is extinguished by dust (e.g. Meurer et al. 
1997). The same appears to be true at high redshift: for ex- 
ample, the brightest few sources in the HDF 850-/im map 
of H98 have very weak or no optical counterparts. If their 
sub-mm emission arises from star formation, all but of or- 
der 1% of the activity is hidden. This is direct evidence that 
the 850-/im background contains a substantial contribution 
from ultraluminous infrared galaxies (ULIRGs) , as has been 
found by a number of other authors (e.g. Barger et al. 1999; 
Lilly et al. 1999). Even where starburst galaxies are visi- 
ble in the UV, measures of the star-formation rate (SFR) 
based on Balmer emission lines often gives rates up to 10 
times larger than would be inferred from the UV continuum 
alone (Pettini et al. 1998). These latter corrections are of- 
ten used to correct the observed UV luminosity density and 
thus estimate the total density of star-formation in the uni- 
verse. The results of H98 show that this is dangerous, since 
there can also be an additional component of star-formation 
that is so heavily embedded that it fails to show up either 
in the UV or in the Balmer lines. The only way to detect 
such a component is through the sub-mm emission, which 
is the approach of this paper. The structure of the paper 
is as follows: section 2 summarizes the optical HDF data 
and shows how to predict sub-mm emission from the HDF 
galaxies as a function of redshift, to within an unknown pro- 
portionality factor; section 3 compares these predictions to 
the observed 850-^im data and shows that there is evidence 
that the 850-pm emission originates over the broad redshift 
range 1 < z < 6; section 4 attempts to test this picture 
by measuring the clustering of the sub-mm emission in the 
HDF map. 
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Figure 1. Comparison of photometric redshifts between the 
Imperial College group (Rowan-Robinson 2000) and Fernan- 
dez-Soto, Lanzetta & Yahil (1999) Solid points show spectroscopic 
redshifts (adopted for the x axis where known). 



2 STARBURSTS IN THE HDF 

2.1 Photometric redshifts 

One great virtue of the HDF data is their wide wave- 
length coverage. This has allowed several groups to esti- 
mate redshifts and spectral types via a variety of template- 
fitting strategies. Generally, these work impressively well 
(see e.g. Hogg et al. 1998 and refs therein). The catalogues 
of Fernandez-Soto, Lanzetta & Yahil (1999; FLY) have been 
particularly successful, as they incorporate JHK data. For 
the present work, we use a set of photometric redshifts con- 
structed according to a similar philosophy at Imperial Col- 
lege (Rowan- Robinson 2000). A comparison of their results 
is shown in Fig. 1; generally, there is excellent agreement. 
The IC database also provides redshift estimates for fainter 
galaxies, beyond the reach of near-IR data. Inevitably, these 
estimates will be less accurate, but are still worth including. 

2.2 UV star- formation rates 

Following Madau et al. (1996), we estimate the star- 
formation rate from the luminosity at a rest-frame wave- 
length of 1500 A: 

SFR L„(1500A) 
Moyr- 1 ~~ lO^-OWHz- 1 ' ^ ' 

A number of slightly different figures for this conversion can 
be found in the literature (e.g. Pettini et al. 1998; Steidel 
et al. 1999); the uncertainty in the appropriate conversion 
factor is however small compared to the uncertain degree of 
extinction. For redshifts 1 < z < 4.4, this wavelength lies 
within the range of the WFPC2 data; for other redshifts, a 
small degree of extrapolation is required. 

The spectral index can also be evaluated (using the two 
bands closest to 1500A): L v oc v~ a . Both the SFR and a re- 
quire the data to be corrected for the effects of intergalactic 
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absorption. We use the corrections given in Fig. 1 of Madau 
ct al. (1996), which may conveniently be approximated by 



= exp [-( 2 /^ c ) 7 ' ] , 



(2) 



Jobs 
/emit 

with z c — (2.50,3.75,4.65) for the (t/300, -B450, Vqog) wave- 
bands. For very high redshifts, these corrections become suf- 
ficiently large that they cannot be made reliably (given the 
uncertainty in redshift). Therefore, for z > 4.5, we assume 
a spectral index of zero and work from the observed magni- 
tude at J814- 

Dust-free starbursts would normally be expected to 
have a ~ 0, whereas many of the HDF galaxies have red- 
der spectra. Pettini et al. (1998) argue that this reddening 
is due to dust. They show that the star-formation rates for 
z ~ 3 starbursts as inferred from Balmer lines are often 
much higher than those inferred from the UV flux, by up to a 
factor 10. Moreover, they claim that the reddest galaxies are 
those for which the correction factor is largest. The spectral 
indices could be used to produce corrected star-formation 
rates, but we shall not do this; the aim of this paper is to 
use the sub-mm emission to infer the total amount of hid- 
den star formation, which can then be compared with the 
activity that is directly visible in the UV. 

Fig. 2 shows the inferred star-formation rate for all HDF 
galaxies with photometric redshifts (assuming £1 = 1 unless 
otherwise stated), plotted against redshift. The uncorrected 
rates are low, with only 5 galaxies exceeding 5 ft _2 M©yr _1 . 
We see that the most active galaxies are found rather uni- 
formly over the redshift range 1.5 < z < 5. To the extent 
that sub-mm flux density is a weak function of redshift at 
high z, this suggests that we would expect a similar redshift 
distribution for 850-^im sources. 

2.3 Dust models and 850-micron emission 

In order to make predictions of the 850-/im flux density, we 
use a simple single-temperature grey-body dust model, with 
T — 50 K and an emissivity proportional to v& with (3 — 1.5. 
The luminosity at rest wavelength Ao is then 



L v 



WHz- 



10 M.o ft -a 



SFR 



y 



-(3+« 



/i^Moyr" 1 J exp(l/y) - 1 ' 



(3) 



where y = Ao/288 /im. The coefficient of proportionality was 
set by comparing with models for Arp220 (Hughes & Dun- 
lop 1999); for a star-formation rate of 75 /i _2 M yr _1 , we 
normalize to an 850- /im flux density of 2.5 mjy at z — 3. As 
is well known, such models predict a sub-mm flux density 
that is very nearly independent of redshift. For z > 1.5, our 
assumptions give an SFR of approximately 30 h~ 2 MQyr~ 1 
for a 1 mjy source. 

Of course, the sub-mm luminosity is a not a measure of 
the total star-formation rate; rather, it measures the fraction 
that is heavily embedded. One way to deal with this would 
be to follow Pettini et al. (1998) and scale the UV SFR by a 
factor that depends on spectral index: SFR — > /(a) SFR; 
the sub-mm emission would then be proportional to the 
corrected star-formation rate, times [1 — l//(a)]. However, 
such a procedure implicitly assumes a simple geometry for 
the dust, and ignores the possibility of an extra contribu- 
tion from heavily embedded star-formation regions. We shall 
therefore prefer to derive the relation between the visible and 




Figure 2. UV SFR against redshift. The raw SFR values 
are low in all cases, with only a handful of objects exceeding 
5fe- 2 M Q yr- 1 . 

hidden SFRs empirically. 

It is clear in advance that very substantial extinctions 
will be required if the HDF starburst galaxies are to be 
detectable in the sub-mm. The rms noise level of the 850- 
/im map is 0.45 mjy, so we require a flux density of 1 mjy for 
an individual detection, or several times 0.1 mjy in order to 
be able to achieve a statistical detection. The former figure 
corresponds to an SFR of about 30 h~ 2 M ( 7 ) yr~ 1 , whereas we 
have seen that the uncorrected UV SFR figures generally lie 
more than a factor of 10 lower than this. In order to measure 
the hidden component, we shall use the term 'predicted' 850- 
fim flux density to denote the value that would be expected 
if the hidden SFR equalled the rate deduced from the UV 
continuum: 



hidden SFR observed 5*850 
visible SFR predicted 5*850 



(4) 



Of course, this makes the assumption that the sub-mm emis- 
sion originates in starlight re-radiated by dust. All that is 
directly measured is the total amount of energy that is being 
re-radiated, independent of its origin. 



3 STATISTICAL DETECTION OF OPTICAL 
STARBURSTS 

3.1 Predicted SCUBA maps 

One direct way in which we can try to see if the general 
starburst population contributes to the SCUBA data is to 
make synthetic maps. Rather than compare directly to the 
observed 850-pim map, we first construct a residual map, 
subtracting the 5 brightest sources (to a limit of 2 mjy) 
studied in H98. The rationale for this is that it is impos- 
sible to predict accurately the brightest few sources, which 
dominate the visual appearance of the map. On the other 
hand, the general background emission is affected by many 
objects, and we stand a better chance of seeing a correla- 
tion between predictions and observations. Fig. 3 shows the 
original and the residual map. 
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Figure 3. The 850-/mi map of the HDF, with North to the top and East to the left. Panel (a) shows the raw map (in which 5 sources 
brighter than 2 mjy have been subtracted and restored with single Gaussian beams, in order to reduce the effects of sidclobcs); panel 
(b) shows the residual map with these 5 sources removed. Both maps are 210" square; the frame has small tick marks at 10" intervals. 
Black denotes high flux; the map is set to zero beyond 100" radius. Both maps have been scaled to remove the trend for instrumental 
noise to increase with radius. 



The residual map (Fig. 3b) shows considerable struc- 
ture. There are a number of isolated regions of peak emis- 
sion around 2 mjy, which are further candidates for point 
sources. Some of these clearly blend, such as the horizontal 
'worm' structure just South of the map centre. This is an 
indication that the image is starting to approach the confu- 
sion limit. There is in addition a hint of a striped pattern, 
with large-scale areas of low and high emission. This pattern 
could be taken as a flat-fielding artefact, and a number of 
reduction strategies were tried, with the conclusion that it 
seemed to be a robust feature. 

We now ask for the predicted appearance of the 850-/xm 
sky, based on the UV starbursts. This is only possible over 
the sub-area of the 850-/im map that overlaps the optical 
HDF, as shown in Figs 4a and 4b. The predicted 850-/im 
map is shown in Fig 4c; this represents a possible future 
view of the sub-mm sky, given an instrument with an ef- 
fective aperture a few times larger than the JCMT. When 
convolved with the JCMT beam (including the effects of 
chopping and nodding), the simulated (noise-free) map of 
Fig. 4d is produced. As a match to the true sky, Fig. 4d 
has some successes and some failures. It does reproduce the 
character of the large-scale variations in emission. This en- 
courages us to believe that features such as the prominent 
'hole' to the bottom right (WF chip 4) do indeed result from 
a general lack of sub-mm sources in that region of the sky. 
Other features match less well, however: the predicted top 
central bright feature is absent in reality; there is no optical 
counterpart to the prominent horizontal blend just South 
of the map centre. This lack of a detailed match is entirely 
to be expected, given that there is a large contribution to 
the map from ULIRGs, which cannot be predicted from the 
HDF UV data. 

3.2 Emission in redshift shells 

It is interesting to dissect the predicted sub-mm emission as 
a function of redshift. It is entirely possible that the relative 



amounts of hidden and visible star formation could evolve as 
a function of redshift, causing one redshift band to dominate 
more than predicted. We therefore divided the map of Fig. 
4d into the contributions from different redshift slices, and 
looked to see how well each correlated with the true map; 
the result is shown in Fig. 5. Although the correlation signals 
from different redshift shells are relatively noisy, there is a 
quite impressive trend with redshift, with little mean cor- 
relation at z < 1, but a positive signal broadly distributed 
over the higher- redshift bins - especially at z > 2.5. Figs 
4e & 4f show the contributions to the predicted sub-mm 
map, divided at this redshift. This is direct evidence that, 
as expected from the SFR data in Fig. 2, the HDF map re- 
ceives contributions from a very wide range of redshifts. The 
UV starburst galaxies thus contribute to the 850-/im back- 
ground in a similar fashion to the ULIRG population, which 
has a median redshift of between 2.5 and 3 (H98; Smail et 
al. 2000). 

Fig. 6 shows in detail the distribution of z > 2.5 galaxies 
(confined to the brighter objects, where the photometric red- 
shifts should be more robust). This is clearly a very strongly 
clustered distribution, as found by Steidel et al. (1998). Per- 
haps the most striking feature is the lack of galaxies in the 
region where the SCUBA map shows its deepest 'hole'. Al- 
though we do not necessarily expect a very precise match be- 
tween the observed and the predicted Ssso for a given object, 
the existence of large-scale features in a clustered galaxy dis- 
tribution will induce structure in the 850-/im map. We will 
return to this point in section 4 below, where we quantify 
the degree of clustering in the sub-mm background. 

3.3 The amount of hidden star formation 

Perhaps the most direct way of seeing if we are detecting the 
effect of the UV starbursts on the 850-/im map is to estimate 
a flux density for each galaxy, and look for a correlation 
between predicted and observed flux densities. We therefore 
take the value of the residual map at the position of each UV 
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Figure 4. (a) The observed SCUBA 850-^tm map, with North to the top and East to the left, windowed over the area of the HDF. The 
frame plotted is 200" square; (b) the SCUBA residual map, with the 5 brightest sources (> 2 mjy) subtracted; (c) the predicted 850-^tm 
map at high resolution; (d) the predicted 850-/im map with the SCUBA beam; (e) and (f) show the same as (d) for respectively z < 2.5 
and 2 > 2.5. Different greyscales are used in each case, to bring out the detail of interest in a given map. The frames show small tick 
marks at 10" intervals. 



source as the best estimate of its true 850-/im flux density. 

These observed flux densities are plotted against red- 
shift in Fig. 7. This is a scatter diagram, with no correlation. 
However, a very different picture emerges when we plot 850- 
fim flux density against UV SFR, as in Fig. 8. Much of this 
plot is also a scatter diagram, but closer inspection reveals 
a clear tendency for the objects of highest star-formation 
rate to have positive fluxes. There is a group of sources with 
Ssso ~ 1 m Jy and SFR > 2fr~ 2 M yr~ 1 that stand out as 



plausible real detections. The parameters for these sources 
are listed in Table 1. Averaging the data in bins yields mean 
values for the flux density at a given SFR, and these results 
are also shown in Fig. 8. This demonstrates that we are 
able to achieve a clear statistical detection of 850-/im emis- 
sion from optically selected galaxies even at a level down to 
0.1 mjy. 

It is interesting to contrast this result with the direct 
approach of pointed observations of individual Lyman- break 
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Figure 5. The cross-correlation coefficient between the 850-^tm 
residual map and the projected UV emission from different red- 
shift shells. The error bars were estimated by rotating the SCUBA 
map relative to the optical HDF data. Although no single redshift 
shell shows a strongly significant cross-correlation signal, there is 
a significant general level of correlation between the SCUBA data 
and the high-rcdshift galaxy distribution, broadly distributed over 
the redshift range z > 1. 
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Figure 6. The brighter (/si4 < 27) starburst galaxies at 
z > 2.5. There is a striking 'hole' in the distribution, centred on 
(x, y) = (40, —40), corresponding to the observed lack of emission 
in the SCUBA map. 

galaxies, which generally yield upper limits of a few mjy 
(e.g. Chapman et al. 1999). The statistical approach clearly 
needs a little thought, since we cannot be certain of the de- 
tection of any individual source. For example, many of the 
features on the HDF map are probably blends of a number 
of sources, so is there not a danger that we might be overes- 
timating the mean flux from any individual source by also 
counting its neighbours? In fact, there is no such bias, at 
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Figure 7. Observed 850-fim flux density against redshift. The 
bottom panel shows the raw values; the top shows mean values 
with standard errors. 
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Figure 8. Observed 850-fim flux density against UV SFR. The 
bottom panel shows the raw values; the top shows mean values 
with standard errors. 

least for a Poisson-distributed field of sources (we consider 
the effect of clustering in Section 4 below). For a Poisson dis- 
tribution, the expectation value of the flux contributed by 
all neighbours is just the mean flux on the map, which is zero 
owing to the chop-nod observing mode. Clearly, the effect of 
blends contributes a confusion-noise term to the general rms 
on the map, and this makes it harder to detect any given 
source. However, the instrumental noise dominates the rms 
in practice, so this is a minor effect. Even though there are a 
finite number of samples of noise on the map, each randomly- 
placed source gives a statistically independent sampling of 
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this noise distribution, so the normal statistics of mean and 
standard error can be applied. 

The mean and standard error of the sub-mm flux den- 
sity to SFR ratio inferred from the points in Fig. 8 is 

SFR 



Ssso/mJy = 0.20 ± 0.04 



/i" 2 M yr- 



(5) 



This implies that the hidden star formation significantly ex- 
ceeds that visible in the UV. A direct fit to the ratios of 
observed and predicted flux densities gives 



hidden SFR 
visible SFR 



= 4.25 ± 1.2. 



(6) 



It is possible that this figure may be an underestimate. 
Examining the entries in Table 1 in more detail, it is appar- 
ent that there are a number of objects with estimated z > 5 
for which (a) the IC and FLY estimates disagree; (b) the 
observed flux densities scatter around zero. The estimated 
redshifts can tend to be unreliable at the extremes, owing 
to detection in only a few bands, so it is possible that the 
sample is corrupted by objects whose apparent high SFR is 
due to an erroneously large estimated z. If the above exer- 
cise is repeated, excluding objects with z cst > 5, the figures 
increase by about la: 



Ssso/mJy = 0.25 ± 0.03 



SFR 



/i" 2 M yr- 



hidden SFR 



6.85 ±1.53. 



(7) 



(8) 



visible SFR 

An alternative approach is to repeat the entire exercise using 
only the FLY redshift estimates; again the effect is to raise 
the level of the statistically detected Ssso by very nearly la. 

It is worth asking if these results are consistent with the 
findings of Chapman et al. (1999), who carried out SCUBA 
photometry of eight Lyman-break galaxies, obtaining only 
one clear detection. The data in their Table 1 gives the ratio 
between 850-pim flux density and SFR, which has a typical 
rms error of 0.2 for an individual object. If fact, the scatter 
of their values is about twice this figure, apparently indicat- 
ing a real scatter in properties. It is therefore appropriate 
to evaluate an unweighted mean and standard error, which 
gives 

SFR 

WmJy = o.i3 ±o.i4 fe _ 2Moyr _ 1 ; 0) 

this is a lower figure than our detection, but is statistically 
perfectly consistent with it. 

Taking the various estimates discussed above, we adopt 
a best guess for the hidden-to-total ratio of 

hidden SFR 



visible SFR 



: 5 ± 1.5. 



(10) 



This is unfortunately a rather model-dependent number, 
since it is sensitive to the spectral shape assumed in the 
sub-mm: the predicted 850- /im flux density scales roughly as 
T~ 3,5 for a fixed bolometric luminosity, so moderate changes 
to the temperature can significantly alter the SFR inferred 
from a given 850- /im measurement. Nevertheless, this correc- 
tion factor for the total SFR is similar to the mean correction 
inferred by Pettini et al. (1998), and other pieces of evidence 
are consistent with their results. Where reliable spectral in- 
dices can be deduced for the sources in Table 1, they are 



large and positive, ranging between 0.3 and 2.0. This is in 
agreement with the claim by Pettini et al. that a red UV 
continuum correlated with a high mean level of extinction. 
Similar arguments were also advanced by Meurer, Heckman 
& Calzetti (1999); they used the observed UV spectra of 
HDF galaxies to predict submillimetre fluxes, and their pre- 
dicted results are in broad agreement with our data. 

However, one should not assume a physical picture in 
which a foreground dust screen extinguishes a dust-free star- 
burst. In a sense, all star formation is dust-enshrouded in 
that stars are born in dusty molecular clouds. Jimenez et al. 
(2000) argue that it takes about 15 Myr for young massive 
stars to burn away this dust and become visible. For contin- 
uous star formation, the hidden-to-visible SFR ratio is then 
really just a ratio of the time spent enshrouded, to the re- 
maining lifetime of the stars once they are visible. Jimenez 
et al. (2000) argue that this ratio should be about 6:1 for a 
Salpeter IMF, and this is not inconsistent with our results. 
It is only possible to exceed this ratio if all star-forming re- 
gions are seen at an early phase, corresponding to a single 
dominant starburst; in this case, the optical counterpart will 
take the form of an ERO (e.g. Smail et al. 1999), or will be 
completely undetectable. 

3.4 Contribution of Lyman-break galaxies to the 
background 

If the true volume density of star formation is about 6 times 
what is visible is the UV, we need to be sure that this is 
consistent with the counts and background at 850 /im. The 
background is related to the proper volume emissivity, j u , 
as a function of redshift by 

(for S7 = 1; see e.g. chapter 3 of Peacock 1999). 

Given the comoving redshift-dependent SFR density, it 
is easy to use this expression to evaluate the background. 
We use the UV density estimated by Steidel et al. (1999), 
with an assumed hidden-to- visible ratio of 5. Integrating to 
z = 5 gives a predicted background of 



z/J® 50 (Lyman - break) = 5.5 x 10~ 10 Wm" 2 sr"\ 



(12) 



which is approximately identical to the total background 
(5±2 x 10~ 10 Wm^sr 1 ) estimated by Fixsen et al. (1998). 
This figure could be increased still further, by a factor of 
about 1.5, if we had adopted the corrections suggested by 
Steidel et al. to allow integration over the whole UV lu- 
minosity function. Alternatively, the estimate might be too 
high by a similar factor, because the SFR density at z < 2 is 
deduced from luminosities at 2800A, and these give system- 
atically larger SFR figures for the more luminous starbursts 
in this paper. Both these corrections are within the statisti- 
cal uncertainty on the SFR density. 

Normal star-forming galaxies thus clearly make a very 
significant contribution to the submillimetre background, 
perhaps the dominant one. At the bright flux densities in- 
vestigated by H98, all sources had hidden-to- visible SFR ra- 
tios that exceeded 60, so ULIRGs dominate the background 
in this flux regime, which contributes roughly 30% of the 
Fixsen et al. background. Active galaxies may contribute up 
to a further 10% (Almaini, Lawrence & Boyle 1999; Gunn & 
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Shanks 1999). Taking the formal uncertainties on the back- 
ground and on the hidden-to-total ratio, the 95% lower con- 
fidence limit on the contribution of UV starburst galaxies 
to the background is 48%, which is just consistent with the 
other contributions. This figure could be lowered in one of 
two ways: (i) if the true background were higher than esti- 
mated by Fixsen et al. (1998); (ii) if this analysis has overes- 
timated the 850-/im flux density that corresponds to a given 
UV SFR. As discussed in section 4, the only way in which the 
second effect might occur is as a result of clustering in the 
background population. In this case, the effective sub-mm 
emission might be overestimated by up to a factor 2. Even in 
this extreme case, however, the UV starbursts must gener- 
ate a minimum of 25% of the 850-/im background. We note 
that this empirical conclusion is in agreement with detailed 
attempts to fit the faint optical galaxy counts with models 
that include a dust correction to the SED, thus predicting 
also the contribution of normal galaxies to the submillimetre 
background (Busswell & Shanks 2000). 

3.5 Comparison with NICMOS data 

To obtain accurate estimated redshifts at z > 4, infrared 
data are required. Part of the HDF (mainly chip WF4) has 
been imaged in JHK by Thompson et al. (1999). On the ba- 
sis of the optical-infrared colours, 8 candidates for galaxies 
at extreme redshifts were identified (Weymann et al. 2000), 
including the z = 5.60 object HDF4-473.0 (NICMOS source 
184.0) studied by Weymann et al. (1998). The distribution 
of these objects is shown in Fig. 9, and their parameters are 
listed in Table 2. Although no significant sub-mm emission is 
associated with HDF4-473.0, there is a striking association 
between the two very close pairs of NICMOS objects and 
the two main 850-/im peaks in the NICMOS field - even 
extending to an agreement in position angle between the 
Northernmost pair and the linear sub-mm source mentioned 
earlier. Given the lack of a plausible candidate identification 
for this object within the optical HDF catalogue, the case 
for association with a z > 5 galaxy seems rather compelling. 

4 BACKGROUND CLUSTERING ANALYSIS 
4.1 Angular power spectrum 

In order to improve on the somewhat anecdotal indications 
of clustering in the 850-/im map, possibly associated with 
high-redshift galaxies, it is necessary to quantify the degree 
of anisotropy seen in the map. The most direct tool for this 
is power-spectrum analysis, which allows a simple decom- 
position of the various contributions to the map. The map 
intensity, /, can be written as 

1 = (S*Wi+N) x W 2 , (13) 

where S is the true sky surface brightness; W\ is the 'window 
function' corresponding to the instrumental beam (includ- 
ing negative sidelobes from chopping & nodding); N is the 
noise, which varies with radius; W2 is a 'censor' function 
that sets the map to zero beyond 100 arcsec radius, and is 
also chosen with a radial variation to make the apparent 
noise independent of radius. 

The structure in the sky background emission can be 
described either by its angular correlation function, w(0), 




Aa / arcsec 



Figure 9. The NICMOS candidate ultra high-redshift galaxies 
are shown as points. Interestingly, these appear to coincide with 
weak SCUBA sources in the case of the two close pairs. 

or by the angular power spectrum, which is best written in 
dimensionless form as Af , the fractional variance per loga- 
rithmic increment of wavenumber: 

Ai J {£6) dl/l 

fao (14) 
Aj=f / w{9) J (ie)9d9. 
Jo 

The angular wavenumber is written as multipole number, 
I, to emphasise that the sky should really be expanded in 
spherical harmonics, Yi m . However, our survey subtends a 
small enough angle that Fourier methods are an excellent 
approximation. For power-law clustering, w(6) = (0/#o) _ % 
this gives 

Ag(l) = (I6 y 2 1 - r( r ( ;/ 2 ) 2) . (15) 

which is equal to 0.77(£#o) e for e = 0.8. This is the canonical 
slope for the small-separation correlations of local galaxies, 
which seems to evolve remarkably little with redshift (Gi- 
avalisco et al. 1998). 

4.2 Window functions 

An interesting aspect of the 850-/im HDF dataset is that 
it contains a number of distinct scales, all of similar mag- 
nitude. The small-scale instrumental noise has an effective 
FWHM of 8.5" the telescope beam has a FWHM of 14.7", 
and chop throws of 30" and 45"; finally, the map radius is 
100", although the noise doubles at radius 90". Each of these 
numbers controls the apparent clustering signal through an 
appropriate window function. 

The window for the beam is most easily evaluated as 
a two-step process. Convolve the sky with a 14.7-arcsec 
FWHM Gaussian, then write the signal as S = S3 — (Si + 
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Table 1 The submillimetre properties of HDF sources with the 
largest star-formation rates (> 2 /i _2 MQyr _1 ). True photometric 
rcdshifts are denoted by an asterisk. 
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u.uy 


1 QA 
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5.81 




2.26 


0.11 


0.85 
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26.89 
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5.29 




2.07 


0.09 


-0.97 


4-439.1 


24.96 


26.03 








4.81 


4.32 


3.07 


0.13 


0.07 


4-445.0 


23.71 


24.08 


24.28 


25.79 


0.50 


2.27* 


2.27* 


2.07 


0.06 


-0.54 


4-454.0 


20.99 


20.79 


20.68 
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0.65 


0.98 




10.00 


0.26 
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0.17 
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2.00 


2.80* 
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2.60 


0.09 


1.37 
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25.22 
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4.67 


4.16 


2.32 


0.10 


0.18 


4-776.0 


24.94 
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5.53 


0.96 


3.80 


0.17 


0.53 


4-878.11 


23.22 


23.53 


23.93 


25.48 


1.10 


2.33 


0.00 


3.24 


0.10 
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Table 2 The submillimetre properties of HDF/NICMOS candi- 
date high-rcdshift galaxies. 



Name 

96.0 
118.0 
123.0 
150.0 
184.0 
248.0 
266.0 
678.0 



#160 

28.81 
27.86 
28.08 
26.39 
26.80 
28.13 
28.76 
27.93 



Juo S$g/mJy 



28.51 
28.65 

27.21 
26.48 
28.10 

27.60 



-0.15 
1.88 
1.95 
-1.16 
-0.19 
0.10 
1.27 
1.35 



5*2 + S4 + S , 5)/4, where (1,2,3,4,5) label the positions in the 
chop/nod sequence. The variance for this sum is easily writ- 
ten down in terms of the correlation function of the smooth 
sky, which allows the window function to be identified: 



\Wi\ 2 = j + jM£6) - J (2£6) - Jo(3£0) 
+ ImU9) + \j (5£6) + lj Q {6te), 

o 4 o 



(16) 



where 6 = 15 arcsec. 

The noise in the published (slightly smoothed) map has 
an rms of 0.45 mjy in the map centre, with a coherence 
length corresponding to filtering with a Gaussian of FWHM 
8.5 arcsec. The spatial variation in the noise rms is divided 
out, so that the censoring function is 



just W2 = W3, where 
W 3 = 2J1(£R)/(£R); 



(18) 



a good approximation to the effect of the radial variation is 
to model W2 by the expression for W3, but decreasing the 
effective value of 7? to 85 arcsec. 

4.3 Shot noise vs true clustering 

The power spectrum of the observed map is a sum of three 
terms: instrumental noise, clustering of the background 
galaxies, and shot noise arising because the number of back- 
ground galaxies is finite. 

If all sources in the map had the same flux density, the 
shot noise would be just 

A s \ ot = J^, (19) 

where N is the surface density of sources. In the practical 
case where the sources have different flux densities, with 
number count dN(S), the effective value of N is 

2 



N eB = 



[fSdN(S)l 



(20) 



Wi = [1 + (r/90 arcsec) 2 ] " 



r < R = 100 arcsec. (17) 



/ S 2 dN(S) 

For Euclidean counts, dN oc S~'^ 2 dS, and the shot noise 
diverges at high flux densities. It therefore makes sense to 
work with the residual map where the brighter point sources 
have been subtracted. This is possible to a limit of 2 mjy in 
the 850-/im HDF map. 

Blain et al. (1999) have shown that the number counts 
to 0.5 mjy are very close to Euclidean in form, with a nor- 
malization of 



Without the radial variation, the £-space window would be N(> 1 mjy) = 7900 ± 3000 deg 



(21) 
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If we assume that Euclidean counts with dN = AS 5 ^ 2 dS 
hold between 2 mjy and some flux density So , then 



<->o 

The background intensity is 



-f 



(22) 



(23) 



(25) 



/ = / SdN(S) = 2ASq 
I s 

The observed background at 850 /xm is 

i// = 5±2x KT 10 Wm^sr -1 = 0.0033 mjy arcsec -2 (24) 

(Fixsen et al. 1998), so the observed count normalization 
(A = 11850) requires So = 0.31 mjy, implying N e g = 
32800 deg~ 2 . 

It is unrealistic to assume abruptly truncated counts, 
and indeed Blain et al. (1999) present evidence that the 
counts have flattened by 0.5 mjy. Consider therefore the 
model 

dN _ AS~ 5/2 

Is ~ i + (sysy-T ' 

which is a more general variant of the form proposed by 
Barger, Cowie & Sanders (1999). We have no direct con- 
straints on 7, but comparison with the faint cm- wavelength 
counts suggests 7=1, and we adopt this value. The back- 
ground is then / = tvASq 1 ^ 2 , and the surface density above 
1 mjy is 2AS ~ 3 ^ 2 [Sl/ 2 —axcta,n(Sf/ 2 )]. The observed counts 
and background then require A = 16500, So = 0.7 mjy, and 
the effective shot-noise source density is N e g — 92600 deg -2 . 
This is probably a more realistic estimate than the previous 
figure; nevertheless, we have altered N e s by a factor of three 
through different models for the faint counts, so this num- 
ber has to be considered somewhat uncertain. Fortunately, 
it can be estimated directly from the data, as shown below. 



4.4 Results 

Fig. 10 shows the observed power spectrum estimated from 
the 850- /im map. The analysis requires a dimensionless fluc- 
tuation field, so the map (after bright-source subtraction) 
was divided by the expected background from sources at 
S < 2 mjy, which is 0.66 of the total, according to the above 
count model, equivalent to a flux density of 0.54 mjy per 
14.7-arcsec beam. The power spectrum was deduced accord- 
ing to the methods described in Feldman, Kaiser & Peacock 
(1994). These authors showed that power estimates are cor- 
related over wavenumber separations of order the reciprocal 
of the survey size, and this effect is noticeable for I < 2 x 10 4 . 

The data are compared with a model which is a sum of 
instrumental noise, clustering, and shot noise. These intrin- 
sic effects are modified in two ways by instrumental effects: 
the small-scale power is reduced by filtering with the tele- 
scope beam (introducing a factor |Wi|), and the large-scale 
power is reduced because the map is constrained to have 
zero mean flux (Peacock & Nicholson 1991); we cannot de- 
tect structures larger than the map size, so the observed 
power must go to zero at zero wavenumber (introducing a 
factor [1 — |M^| 2 ]). The overall model is thus 

A 2 bs = (A 2 lus + A 2 hot ) x | m| 2 x (1 - |VF 2 | 2 ) 




angular wavenumber SL 

Figure 10. The two-dimensional power spectrum of the 850-^tm 
HDF map. Angular wavenumber is denoted by £, and A 2 (i) is 
the contribution to the fractional sky variance from unit range 
of ln£. The dotted line shows the expected noise spectrum; the 
dot-dashed line is the shot-noise component for an effective sur- 
face density of N e g = 92600 deg -2 (see text); the dashed line is 
the contribution from clustering with Bo = 1 arcsec; the solid line 
is the total model power, which fits the data well. 

Shot noise is distinguishable from true clustering be- 
cause these effects have a different dependence on scale. 
Our data cover only a restricted range of scales, so in prac- 
tice we shall assume that the true angular clustering has 
a correlation structure of the same shape as is observed for 
most faint galaxy samples, but with an unknown amplitude: 

w(0) = (e/e () r - 8 . 

Fig. 10 shows an illustrative fit of shot noise with N c g = 
92600 deg~ 2 as estimated above, and clustering with 0q = 
1 arcsec. However, both these parameters can be fitted to the 
data in the usual maximum likelihood manner. The result of 
this exercise is shown in Fig. 11. As would be expected, there 
is a degree of degeneracy in the free parameters; the data 
indicate an excess of power in the residual map above pure 
noise for 10 4 < I < 5 x 10 4 , and this constrains the sum of 
clustering and shot power. The different scale dependencies 
of the two signals mean that the degeneracy is broken, and 
the preferred model has significant clustering (Go — 3 arcsec) 
and low shot noise (N e ff — 5 x 10 5 deg -2 ). However, this is 
an unrealistically high surface density; we should really fix 
N e g at its a priori estimated value, in which case there is no 
significant detection of clustering, and 9o = 3 arcsec can be 
ruled out at approximately 95% confidence. 

4.5 Implications for origin of the background 

What is the interpretation of the relatively weak angular 
clustering? The angular power is a projection of the true 
3D spatial power spectrum, given by the k-space version of 
Limber's equation: 



Al 



-if 



& 2 Wy) C(y)y 5 <f> 2 (y) dy 



(27) 



+ A 



x(l-|^ 3 | 2 ) 



(26) 



(Kaiser 1992; see section 16.6 of Peacock 1999). In this ex- 
pression, y is comoving angular-diameter distance, and the 
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2x10' 



5x10" 10" 2xlO u 
N e „ / deg- 2 



5xlO J 



Figure 11. Contours of likelihood relative to maximum for the fit 
of the 850- /im power spectrum data by a combination of shot noise 
(effective surface density for sources with flux density < 2 mjy 
iV e ff) and galaxy clustering (where 8o is the scale-length in the 
angular correlation function). The contour interval is 0.5 in ln£, 
corresponding to confidence contours of approximately 39%, 63%, 
78%, 87%, 92% etc. 



'selection function' cj>(y) is normalized so that 

poo 

/ y 2 4>(y)C{y)dy = 1. 
Jo 

The curvature factor 



C(y) = 



1 + 



(i - ft)y 

c/H 



-1/2 



(28) 



(29) 



so it is unity for spatially flat models. 

If A 2 (k) oc fe 1 ' 8 (so that A 2 e (£) oc <? ' 8 ), the angular 
power scales as 

Ag oc hl^Ci^ 
[jy 2 ^y)C(y) dy] 

The redshift distribution relates to the selection function via 
n(z)dz oc y 2 (f>(y)C(y)dy. What this says is that the ampli- 
tude of angular clustering depends on the range of redshifts 
that contribute to the sample under study, and a low value 
of 8 requires a large range of redshifts to be contributing. 

To be quantitative, consider the results of Giavalisco 
et al. (1998) for the angular clustering of Lyman-break 
galaxies. They deduce 8 = 2 arcsec over a redshift range 
2-5 < z < 3.5. The relation of comoving distance to observ- 
ables is 



C(y) dy=—dz x 
Ho 



[(1 -Q.ru- 0,)(1 + zf + Q V + Q m (l + zf] 



-1/2 



(31) 



and this allows us to see what would happen to the angular 
clustering signal if the redshift range is extended. Assuming 
a uniform contribution to the sub-mm surface brightness 
from 2 < z < 6, as indicated by the results in earlier sections, 
the angular clustering signal would be expected to drop by 
a factor of 0.35 (fi m = 1, Q v = 0) or 0.34 (fi m = 0.3, 
f2„ = 0.7) compared to that measured by Giavalisco et al., 



for the same degree of spatial clustering. 

In making this comparison, it is important to be clear 
that we are including the angular clustering of the ULIRG 
population, rather than just that of Lyman-break galaxies, 
since the former give a significant contribution to the back- 
ground. Nevertheless, the Lyman-break galaxies serve as a 
useful reference. They are a highly biased population, as in- 
deed would be expected for all reasonably massive galaxies 
at high redshift (Steidel et al. 1998; Peacock et al. 1998). 
It would therefore not be surprising if the ULIRG popu- 
lation displayed a similar level of spatial clustering. How- 
ever, even if this is so, the sub-mm background arises over 
a broad range of redshifts, and the sky at 850-/im is there- 
fore expected to be more uniform than the Lyman-break 
sky. Strong clustering of the 850-/im background would only 
arise if the ULIRGs were much more strongly clustered that 
Lyman-break galaxies; our results indicate that this cannot 
be the case. 

As a final point, we can now address the question of 
whether clustering in the background could bias the statis- 
tical detection of flux from the UV starburst galaxies. Any 
given galaxy will statistically be surrounded by a 'halo' of 
emission contributed by correlated galaxies, whose surface 
brightness at angular separation 6 will be w(8)I, where I is 
the mean background. Since we are dealing with two pop- 
ulations, UV starbursts and ULIRGs, in principle we need 
both the UV-ULIRG cross-correlation function as well as 
the UV autocorrelation function. However, since a correla- 
tion coefficient must be less than unity, the cross-correlation 
cannot exceed the larger of the two autocorrelations. We can 
therefore get a limit to the effect by assuming 8o = 2 arc- 
sec. The effective point-source flux is obtained by integrating 
w(8)I v over the beam. The result of this exercise is an ef- 
fective point-source flux density of 0.19 mjy. This is large 
enough to be a dominant bias for the starbursts with SFR 
around 1 /i _2 MQyr" 1 , but cannot erase the signal from the 
more luminous objects. In order to obtain an upper limit to 
the overall bias that could be introduced by this effect, the 
analysis of section 3.3 for the mean level of 850- /im emis- 
sion was repeated, subtracting 0.19 mjy from all the flux 
densities. The result is to reduce the mean sub-mm output 
for a given SFR by a factor of approximately 2, although 
we have argued that this is very much an upper limit to 
the correction. The effect of clustering can therefore in prin- 
ciple have a significant quantitative impact on issues such 
as the total sub-mm emission from UV starbursts (section 
3.4). However, this does not alter the main point, which is 
that UV starbursts are one of the largest contributors to the 
background radiation at these frequencies. 

5 CONCLUSIONS 

In this paper, we have investigated whether the SCUBA 
850-/im map of the HDF can set useful constraints on the 
general population of optically-selected starburst galaxies, 
or whether we are confined to investigating only totally ob- 
scured ultraluminous infrared galaxies. Although data of im- 
proved resolution and sensitivity will be required in order to 
study the sub-mm properties of HDF galaxies on an individ- 
ual basis, we believe that there is good statistical evidence 
that the galaxies detected in the optical HDF catalogues are 
emitting in the sub-mm band, with a flux density of about 
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0.2 mjy for an apparent UV star-formation rate of 



This level of emission is consistent with the 



Ivison R.J., Smail I. 
A.W., Owen F.N. 



idea hat the UV emission from these galaxies underesti- 



mates the total star-formation rate in most cases, by a mean 
factor of approximately 6. This means that the Lyman-break 
population must contribute at least 25% of the background 
radiation at 850 /im. 

We have used a variety of tests, including cross- 
correlation with UV data and power-spectrum analysis of 
the 850-pm data, to demonstrate that the sub-mm back- 
ground receives approximately equal contributions from a 
very wide range of redshifts, 1 < z < 6. Indeed, some of 
the strongest individual examples of evidence for an associ- 
ation between high-z galaxies and features in the sub-mm 
map involve galaxies towards the higher end of this range. 
Together with the lack of a clear optical counterpart for the 
brightest HDF source (HDF 850.1; Downes et al. 1999), this 
suggests that the upper bound of the redshift distribution 
for 850-/im sources may not yet have been reached. 
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